An integrated technology for quantitative wide mutational scanning of human antibody Fab libraries

Antibodies are engineerable quantities in medicine. Learning antibody molecular recognition would enable the in silico design of high affinity binders against nearly any proteinaceous surface. Yet, publicly available experiment antibody sequence-binding datasets may not contain the mutagenic, antigenic, or antibody sequence diversity necessary for deep learning approaches to capture molecular recognition. In part, this is because limited experimental platforms exist for assessing quantitative and simultaneous sequence-function relationships for multiple antibodies. Here we present MAGMA-seq, an integrated technology that combines multiple antigens and multiple antibodies and determines quantitative biophysical parameters using deep sequencing. We demonstrate MAGMA-seq on two pooled libraries comprising mutants of nine different human antibodies spanning light chain gene usage, CDR H3 length, and antigenic targets. We demonstrate the comprehensive mapping of potential antibody development pathways, sequence-binding relationships for multiple antibodies simultaneously, and identification of paratope sequence determinants for binding recognition for broadly neutralizing antibodies (bnAbs). MAGMA-seq enables rapid and scalable antibody engineering of multiple lead candidates because it can measure binding for mutants of many given parental antibodies in a single experiment.


Supplementary Figures
Supplementary Figure 1 | Schematic of assembly strategy with shuttle vector and yeast display mutations.a. Yeast display plasmid map highlighting most of the relevant genes (shown is kappa only; the lambda map is otherwise identical).The gene segments are not drawn to scale.b.Two-step cloning strategy for assembling barcoded Fab libraries.Along with a bidirectional promoter (BDP) plasmid, any VH-VL pair or library is assembled by Golden Gate into a minimal 3.6 kB mutagenesis plasmid containing a CamR selection marker, a high copy number ORI, and regions of homology to the CH1 and light chain sequence.There are small mutagenesis destination plasmids for both human kappa and lambda light chains.After mutagenesis, antibody mutant libraries are assembled by the method of Gibson with a yeast surface display plasmid containing a unique UMI.This final plasmid is identical to that in panel a. c.Summary of mutations in the yeast display vector.The major missense change is S5Q on CH1 necessary for encoding a PaqCI restriction site near the CDR H3 for short read sequencing pairing of UMI to the VH gene.We also removed the light chain 6x Histidine epitope tag and replaced it with a V5 epitope tag (GKPIPNPLLGLDST) to be able to measure binding to antigens that may be His-tagged with an anti-His conjugated fluorophore.Unique PaqCI and BbsI sites are necessary for UMI-Fab pairing by short read sequencing; antibody sequences encoding these sites are not compatible with the short read sequencing protocol.Three isogenic plasmids (one CC12.1 variant, two 4A8 variants) were mixed at different molar ratios and the intramolecular ligation for the VL chain and amplicon protocols were performed in triplicate (n=3).We varied the following parameters: polymerase, number of PCR cycles, and PCR clean-up method (rSAP (New England Biolabs) or column cleanup).Amplicons were deep sequenced on the same MiSeq flow cell, and observed frequencies of pairing were extracted.True frequencies of correct pairing (freqtrue) were obtained for the lowest abundant variant using the following equation:

Supplementary
Here,  )*+ is the observed frequency by deep sequencing and f is the fraction of the lowest abundant variant.As this fraction goes to zero the true frequency is identical to the observed frequency.True frequencies and p-values from paired, one-tailed t-tests are reported.We also report true frequencies and p-values of performing the intramolecular ligation separately on individual plasmids (IML control).The protocol chosen for barcode-Fab haplotyping is highlighted in green.

. INFERRING ANTIBODY-ANTIGEN BIOPHYSICAL PARAMETERS FROM SEQUENCING DATA: EXECUTIVE SUMMARY
In deep mutational scanning, a population of mutational variants of a protein is passed through a selection or screen; this screen changes the underlying frequencies of each of these variants.Deep sequencing is used to count each variant in the population, which is used to infer the frequency of each variant in the population in a reference population and after the screen.This frequency change is converted to some score that, ideally, relates to the functional properties of the variant.This technical note describes our framework for inferring, from this processed sequence data, both dissociation constants and maximum fluorescence for antibody variants encoded in yeast displayed protein libraries screened by fluorescence activated cell sorting (FACS).
In FACS, populations are screened by collecting cells with fluorescence above a certain gating fluorescent threshold, or between two fluorescent gating thresholds; cells sorted according to these fluorescent gates are said to be sorted into bins.A clonal population of cells will exhibit a mean fluorescence with a certain variance according to cell size, surface density of displayed proteins, or other factors.Thus, only a fraction of cells for each variant will exceed the fluorescence threshold needed for collection into a given bin.Furthermore, if the fraction of cells that are sorted into a bin is known, one can infer the likely mean fluorescence for a given variant at that labeling concentration.Finally, as described in further detail below, sequencing data and other experimental observables can be used to infer the fraction of cells collected by the gating strategy and thus the mean fluorescence of a variant for a given labeling concentration.Some of the descriptions below come in part from Kowalsky et al 1 , and Kowalsky et al 2 .
We seek to infer variant-specific dissociation constants (Kd,i) using, for example, the Hill equation below: Here is the mean fluorescence of cells displaying variant i at a given labeling concentration [L], Fmax,i is the maximum fluorescence for the variant i, and Fmin is cellular autofluorescence.
Intuitively, if we can infer the mean fluorescence at different labeling concentrations ([L]), we can reconstruct isothermal titrations for each variant i (e.g. , vs. [L]) to find a best fit Kd,i and Fmax,i using non-linear regression.An example from barcode ATGCACACATTTAAAGCTGT corresponding to variant 4A8 M59I is shown below in Fig Note S1.
We can approach this inference problem by regression, as it allows for the quality of the model fit to the data using the chi squared metric while also giving robust methods for confidence interval testing.As will be shown, we can also use maximum likelihood estimation in a quantitatively identical way.However, we cannot regress on the reconstructed mean fluorescence, as error is not distributed uniformly in both directions.Instead, we regress on the vector of probabilities of sorting into each bin.
Supplementary Figure 11: Fluorescence reconstruction for barcode ATGCACACATTTAAAGCTGT at 10 labeling concentrations.
In summary, we label our population of antibody variants at different antigen concentrations and use FACS to sort these antibody variants into different bins.Sequencing of these populations allows us to reconstruct the likely mean fluorescence for a given labeling concentration by inferring the fraction of each variant that is present in a binned population.Summing over all labeling concentrations allows us to find the most likely parameter value for dissociation constants, the confidence interval associated with that parameter estimation, and the quality of the fit using weighted nonlinear regression.

WHAT IS THE PROBABILITY OF A CELL COLLECTED ABOVE A CERTAIN FLUORESCENT THRESHOLD?
Let's call this fluorescence threshold a gating fluorescence (Fg) and ask for the probability that a given clone i exhibiting a mean fluorescence intensity ( , ) will be captured by this gate.A graph of this relationship is below: :√7 / (7)

WHAT DO THESE PROBABILITIES LOOK LIKE IN PRACTICE?
For a typical monovalent binding experiment, one labels yeast cells displaying a binding protein with a fluorescently conjugated ligand.We have found that  for phycoerythrin (SAPE) labeled populations range from 0.9-1.05 2 .Let's assume a  = 1.02 for this example.We find that many protein-ligand interactions we consider in lab to be well fit by a Hill equation with no cooperativity: For the phycoerythrin (SAPE) labeled populations we usually consider, a typical value of  -,1 representing cell autofluorescence is 350 MFI in our experimental set-up using a Sony SH800 cell sorted with a 488 nm laser and compensation for fluorescein.The two protein-specific terms are the max fluorescence (Fmax,i) and the dissociation constant Kd,i for the interaction.These will be variant-specific.For reasonably expressed and well-behaved proteins our Fmax,i is typically in the 50,000 MFI range.Let's nondimensionalize the ligand concentration so we can remove one variable.

WE CAN INFER THE PROBABILITY USING FREQUENCY DATA
Our observable for deep sequencing experiments is a set of read counts for variant i in each bin j and for each labeling concentration k (let's call these read counts rijk).Additionally, we have the reference read counts we can observe for variant i (let's call this rir).We can directly convert observables to probabilities of sorting into a given bin j by comparing these read counts to those from the reference population.The reference population is critically important given that the comparison is the probability of being captured by a gate relative to the condition of no gate.Therefore, your reference population must be identical except for the fluorescence gate you sort at.
We write the probability as the number of cells of variant i collected in the j th bin and k th labeling concentration (xijk) relative to the number of cells of variant i in the reference population (xir): The frequency of variant i (fijk, fir) is just the number of counts observed divided by all counts, so we can write: Here ∅ is the total fraction of cells collected in the sorting bin relative to the reference population, and the frequency of each variant has been converted to experimental observables derived from deep sequencing.Equation (10) is the second fundamental equation because it states that the probability pijk (set by Fmax,i and Kd,i) we observe for a given labeling concentration k and bin j are a function of the observables from the deep sequencing experiment.

SOURCES OF NOISE IN RECONSTRUCTING FLUORESCENCE FROM EXPERIMENTAL OBSERVABLES
A major challenge for sequence-function reconstruction experiments comes from determining the appropriate confidence level set for each experimental measurement.This is important as low and high values of pijk give large uncertainties in the measurement of  ,= (see Fig Note S3 below).
Supplementary Figure 13: Mean fluorescence of variant as a function of pijk -the probability of sorting into a bin above some Fg = 30,000 MFI.Low and high observed probabilities result in large, one-tailed uncertainties in the value of the mean fluorescence.
For parameter inference it is important to identify and quantify sources of noise in the fluorescence reconstruction.Intrinsic noise comes in the act of sorting discrete cells, preparing amplicons from yeast cells by PCR, and sequencing discrete nucleic sequences.Extrinsic noise results from the efficiency of the overall process of cell sorting and recovery.
We have previously shown 5 that deep sequencing read counts from FACS data can be evaluated according to Poisson probability distributions.We have previously determined the propagation of errors for Poisson noise in the counts of the reference and selected populations 6 .Although this is an underestimate of error because population bottlenecks occur during cell sorting, outgrowth, and amplicon prep leading to overdispersion, we find empirically that Poisson noise is a reasonable approximation for well-designed experiments.By propagation of errors, we can determine the variance associated with the value of  ,<= : >,<=,,1!",1+,?

7
=  ,<= (∅ >,<=,,1!",1+,? The other source of error is extrinsic relating to error rate in the sorting itself -what is the probability of a mis-sorting event?This appears to be 2% in our experimental set-up, but we expect that this error rate may vary. >,<=,$/!",1+,?7 = (0.02) 7 (15) In the experiments presented in Figures 4 and 5 of the main text, we included Fab nonbinders to measure this error directly from the sequencing data.For these experiments, these values were observed to be 0.014, close to the values used in the initial experiment.

PARAMETER ESTIMATION USING MAXIMUM LIKELIHOOD ESTIMATION
The log likelihood framework states that the parameter set most likely to fit a given set of data occurs with maximization of the summation of the log probabilities of each experimental measurement: , ( A,, ,  -./,, ) = (∑ ' ,<= ) Here, Modelijk is the model probability (given parameters  A,, ,  -./,, ) of a variant i being sorted into bin j at labeling concentration k.We must assume some probability distribution -given the sources of noise and the fact that reference and sorted counts are typically >10, a Gaussian probability distribution is justifiable here.Expanding terms, we can write: , ( A,, ,  -./,, ) = (∑ ' ,<= ) Expanding the Gaussian probability distribution and removing constant terms, we arrive at: , ( A,, ,  -./,, ) = ∑ − (20)

CONFIDENCE INTERVALS USING MAXIMUM LIKELIHOOD ESTIMATION
Using a MLE framework that minimizes the chi squared metric ( -,1 7 ) results in the simplification of confidence interval measurements.As such, we follow standard approaches 7,8 for determining 95% confidence intervals using the critical value of the F distribution statistic (F0.05) using the following equation: Where n is the number of experimental data points (here, the number of bins used for MLE), and  7 is the chi squared metric for given parameter values of Kdi and Fmaxi.

Figure 2 |
Yeast surface display titrations.Isogenic yeast surface titrations for antibodies reported in main text Figure 1c.Error bars represent 1 s.d. of 2 measurements.The curve fit shown is a Hill equation where the Hill coefficient is constrained to unity.Supplementary Figure 3 | Intramolecular Ligation Products.Intramolecular ligation (IML) reaction products from 1 µg of barcoded 4A8/CC12.1/COV2-2489Fab library (lanes 6-9) were separated by a 1% (w/v) agarose gel.Lanes 2-4 show individual controls reactions without ligase.The ligation product from BbsI is 1.8 kB, while that from PaqCI is 6.4 kB.In these reactions, the UMI is paired to the VL with BbsI intramolecular ligation and to the VH with PaqCI.Biological replicates (Rep1, Rep2) were performed, yielding reproducible band sizes and intensities.The gel shown is the unmodified and uncropped image.Supplementary Figure 4 | Optimization of PCR amplicon preparation for barcode-Fab haplotyping.

Figure 5 |Supplementary Figure 6 |
Repeated 'CGGCGG' motif in COV2-2489 sequences causes drop in quality at nucleotide 147 on VH reverse read.Average quality score vs. sequence position for 4A8 & CC12.1 antibodies (orange) compared to COV2-2489 (blue).The inset shows the nucleotide sequence adjacent to the drop in quality score.Cytograms from first sort with S1 and 4A8/CC12.1/COV2-2489Antibodies.Cytograms showing sorting gates for first demonstration of the method with mixed Abs against S1.Cells were first gated for yeast cells, single cells, and cells displaying Fabs before being gated and sorted for the top 25% and next 25% bins.Supplementary Figure 7 | MLE KD estimation for grouped barcodes with 95% confidence intervals.(Top row) Histograms of MLE KD estimates for each barcode with calculated mean absolute error from isogenic titration data.(Bottom row) 95% confidence intervals for each barcode (blue X, S7T: n=30, M59I: n=12) and from barcodes collapsed by variant (orange X)

Figure 10 |
Potential development trajectory for SARS-CoV-2 antibody 2-7.(a) Sampling of 6 of the 30 potential intermediates between the UCA and mature 2-7.The affinity of each variant is shown as ΔΔG (kcal/mol) relative to the mature 2-7 antibody.Mature 2-7 has an inferred Kd of 9.6 nM and the UCA a Kd of 255 nM.(b) LASSO regression one body weights for ΔΔG for the five VL mutations.Weights in kcal/mol are shown relative to the mature 2-7 Ab.The E52D mutation is energetically unfavorable and unlikely to have appeared except in conjunction with the N55K mutation.

:
this expression is equivalent to minimizing the weighted sum of square errors or the chi squared metric.The algorithm changes the probabilities of Modelijk by changing parameters in the Hill function, and we use off-the-shelf optimization software to find the minimization of the function.−, ( A,, ,  -./,, ) = ∑ (> &*+ 'C)A$8 &*+ E.EG ( − 1, )) (21) , , ) =

Supplementary Table 1 reports
the resulting probability lookup table:This table shows that the large standard deviation resulting from square gating gives useful probabilities at many different gating fluorescence values representing different dissociation constants and/or max fluorescence values.